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Abstract 

We report our results for the bag-parameter Bk obtained from the quenched sim- 
ulations on the lattice with Wilson fermions at three values of the lattice spacing. 
We implemented the method by which no subtraction of the mixing with other four- 
fermion AS" = 2 operators is needed. Our final result, in terms of the renormalisation 
group invariant bag-parameter, is Bk = 0.96 ± 0.10. 



PACS: ll.15.Ha, 12.38.Gc, 13.25.Hw, 13.25.Jx, 13.75Lb, 14.40.-n. 



1 Introduction 



The parameter characterising the size of the indirect CP violation in the system of neutral 
kaons, Ek, has been accurately measured since long ago 0. The precise theoretical estimate 
of the corresponding K°—K° mixing amplitude, however, is still missing due to uncertainties 
in the computation of the matrix element of the operator 

Q AS=2 = ( -A^ (1 _ 75 )^ )( -^ (1 _ =Ql+Qu (1) 

where Qi and Qi are respectively the parity conserving and parity violating part of AS=2 . 
A and B are the color indices. The matrix element of the renormalized operator, 

(K°\6 AS =^)\K°) = (K'lQ^lK ) = \fW K BM , (2) 

is conveniently parameterised in terms of the bag-parameter B^, the measure of the devia- 
tion of the matrix element from its value obtained in the vacuum saturation approximation 
(in which Bk — 1). 

Over the past two decades quite an impressive progress in computing Bk on the lattice 
has been made. We now know how to renormalize the four-fermion operator Qi non- 
perturbatively in the RI/MOM [2] and in the Schrodinger functional scheme jSj- We also 
know how to relate to Bk, the renormalisation scheme invariant bag parameter, since the 
anomalous dimension is calculated in a number of schemes at NLO in continuum perturba- 
tion theory [U E| , the same accuracy at which the corresponding Wilson coefficient has 
been calculated [7.. A high statistics computation of Bx with Wilson fermions for several 
values of the lattice spacing a has been performed in ref. jH]. The preliminary unquenched 
calculations have been made too [Oj. However, all the works in which the Wilson quarks 
were used suffer from the potentially large systematic uncertainty that arise from the large 
mixing of the operator Q\ with other parity-even operators Q2-5 which have different naive 
chiralities. 1 That feature is a consequence of the explicit chiral simmetry breaking induced 
by the Wilson term in the quark action. In other words, the renormalization pattern of the 
lattice operator Qi regularized a la Wilson is 

5 

= Z(a/J,) Qi(a) + y^ j A i (a)Q i (a) , (3) 

i=2 

where Z(afi) is the multiplicative renormalization constant present also in formulations 
where chiral simmetry is preserved, while A 2 _5(a) are the mixing coefficients peculiar for the 
Wilson regularization. The difficulty is not only that one needs to compute the subtraction 
constants A 2 -5(a) non-perturbatively but one should also compute them very accurately 
because the lattice regularized bare matrix elements (Q2-5) are orders of magnitude larger 
than (Qi). Therefore, even though the subtraction constants are numerically very small 
(see the tables in ref. [HI| 2 ), the net effect of the subtractions is large. It is clearly desirable 

1 For the explicit forms of all the parity even operators, Q1-5, see e.g. Section 3.2 of ref. 
2 To avoid a notational ambiguity we point out that the subtraction constants Ai2(a), . . . , Ais(a) pre- 
sented in ref. ^D] correspond to ^2(0), . . . , ^5(0,) used in this letter. 
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to have a method that allows one to compute the matrix element (J2J) without necessity to 
subtract the mixing. 

In this letter we use the hadronic Ward identity, proposed in ref. ^l] , to relate the matrix 
element of the operator Qi to the parity violating one, Q\. The latter does not suffer from 
the spurious mixing and thus the problem of mixing with other dimension-six operators is 
circumvented. The price to pay is that one has to compute a four-point correlation function 
where one pion is integrated over all lattice space-time coordinates. Similar in spirit, but 
quite different in practice, is the proposal made in ref. J2| where the chiral rotation has 
been added to the mass term as to kill out the spurious lattice mixings. Preliminary study 
of the Bk parameter by using that method, has been presented in ref. ^3] ■ 

In sec. |21 we will briefly recall the basic elements of the Ward identity method to compute 
Bk without subtractions; in sec. El we present the results for the matrix element (J2J) for the 
directly accessible pseudoscalar meson masses from which we will extract the B^ parameter; 
in sec. EJwe briefly conclude. 



2 Strategy 

In this section we will briefly recall the main steps involved in the extraction of the Bk- 
parameter by using both methods, with and without subtractions. 



2.1 Standard Procedure: "with subtractions" 

The standard way to extract the matrix element (J2J) proceeds through the computation of 
the correlation functions 

G K o(t) = (K° p \t)K P (0)) , %(4,y = {K%\t x )QMKl{t y )) , (4) 

with Qi defined in eq. ©, and K P (t x ) = ^2 x -d(x) , y5s(x). Therefore to get GU (t x , t y ) one 
must compute the correlators by using the complete set of parity conserving four fermion 
operators, Qi-$, subtract the spurious mixing, and provide the multiplicative renormalisa- 
tion, as indicated in eq. El This procedure is particularly delicate because the approximate 
restoration of chiral symmetry (which is exactly recovered only in the continuum limit) 
depends on how well the subtractions are made. The subtraction constants A 2 _5 do not 
depend on the renormalization scheme. Their values have been recently estimated non- 
perturbatively, in the RI/MOM scheme |lUj . 

The matrix element (J2J) is extracted from the study of the large time asymptotic be- 
haviour of the ratio 

™w rn G Qi (t x ,t y ) (^°|Qi|if°) () 

{y) " G K o(t x )G K o(t y ) *■»*■»*/» \(0\K p \K°}\*> {> 

where we fix one of the source operators at t x so that the kaon state which is created by the 
four fermion operator in the origin is already asymptotic when annihilated by K p (t x ). The 
time t y , instead, is free. On the plateaus, T ^> t y ^> T/2, where all the operators are far 
away from one other, we read off the desired matrix element divided by the pseudoscalar 
density squared. 
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2.2 Alternative Procedure: "without subtractions" 

The method proposed in ref. ^l] is based on the use of a Ward identity which arise from 
applying the T3 axial rotation, 

Su(x) = ia{x)^^u{x) , 5u(x) = ia{x)u{x)^^ , 
5d(x) = — ia(x)j5d(x) , Sd(x) = —ia(x)d(x) r y 5 , (6) 

onto the matrix element (K P (x)Qi(0)K P (y)) , where K P = ZpK P . To write down the 
relevant Ward identity, we introduce the bilinear operators 

K° s (t) = d(x)s(x) , n°(x) = d(x)-f 5 d(x) - u{x)^u(x) , (7) 

X 

and the corresponding renormalized Kg(t) = ZsKg(t). With these definitions in hands the 
renormalized Ward identity reads 

2(k P (t x )Q 1 (0)k° P (t y )) = 2mJ2(^ O ^)K°p(t,)Qi(0)k o P (ty)) 

z 

-(k s (QQi(o)kUty)) - (^(^)4(0)^(^)) + 0(a) , (8) 

where, in view of the fact that we work out of the chiral limit, we dropped the sum over 
the space-time of the term containing the divergence of the axial current, which appears, 
together with the first term on the r.h.s. containing 2mU°(z). The dropped term is zero 
when there is no momentum transfer between the initial and final states, provided that no 
singularities are encountered. Working at non-zero quark mass then ensures that infrared 
divergences are avoided, so that the integration over space-time indeed yields a vanishing 
value. In practice, we work in the SU(3) limit, i.e., we take all three quarks to be degenerate 
in mass, m u = irid = m s = m. The term on the l.h.s. of eq. (jHJ), corresponding to the rotation 
of the operator Qi, is the desired matrix element. The last two terms in eq. (jHJ) correspond to 
the rotation of the pseudoscalar kaon sources. These terms, although necessary to saturate 
the Ward identity, disappear in the SU(3) limit as shown in appendix. Thus, the Ward 
identity we use in practice reads 

CK£(t x )Qi(0)#£(t tf )) = m(afi)Z(afi) ^(fl ^)^ (t x )Q 1 (0)K° P (t y )) = G Ql (t x ,t y ) , (9) 

z 

where we stress the presence of 0(a) artefacts, i.e., the four fermion operators are not im- 
proved. Owing to CPS symmetry the parity-odd operator Qi renormalizes multiplicatively 
only. Z{a[i) has been recently computed non-perturbatively in the RI/MOM scheme |10j . 
We use the quark mass, m(a/i) = pZ J 4(a)/Zp(a/x), defined through the non-singlet axial 
Ward identity, 

<W)*j?(0)> 
2<Kj(t)tf?(0)) 

where A^it) = ^ x s(x)7 M 7 5 c?(x), and Za(o) is the axial-current renormalization factor [10J. 
Notice also that in eq. (jHJ) the renormalisation constant of the pseudoscalar density, Z P (afi), 
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in m(afi) cancels against the one in ft (z). In terms of Feynman diagrams, eq. Q can be 
written as 

CS{t x ,t y )+D8{t x ,t y )j = 2Z A p\cE{t x ,t y ) + CE(t y ,t x ) + DE(t x ,t y ) + DE{t y ,t x )j , (11) 

where C8(t x ,t y ) and D8(t x ,t y ) correspond to the connected and disconnected "eight" dia- 
grams, while CE(t y ,t x ) and DE(t y ,t x ) refer to the connected and disconnected "emission" 
diagrams shown in fig. 1 of ref. jllj . Proceeding like in the standard method, the matrix 
element is extracted from the study of the ratio, 



^j>w/o subtr. ^j. 



(K°\Qi\K°) 

G K o(t x )G K o(t y ) \(0\K° P \K°)\* ' 



(12) 



3 Extraction of Bk 

In this section we present our main results. We use both procedures, the standard and the 
one without subtractions, which provides us a useful cross-check. Of course the two methods 
suffer from 0(a) effects that are different in size, but should converge to the same value in 
the continuum limit. 



3.1 Lattices and signals for the ratios (03) and (1121) 

The details of our lattice setups were presented in our previous publications (TU1 E]. We 
work at three lattice spacings which correspond to /3 — 6.0, 6.2, and to 6.4. In each 
simulation we work with four different values of the quark mass, i.e., with four values of the 
parameter k, and compute the correlation functions needed to form the ratios © and ()12j) . 
In fig. ^ we show the signals we obtain by using both methods and for all quark masses used 
in our simulations (/c's are ordered as mi > m<i > > m^). The plateaus correspond to 
the signals for the bare operators, i.e. without accounting for the overall (scale dependent) 
renormalisation constants Z{a\x) and Z{a\x). For the standard method we need to specify 
the subtraction constants A 2 _5(a). We use the results recently obtained in ref. [10J. In our 
calculation one source operator is fixed at 



t x = 12 


,15 


,17 






0=6.0 


13=6.2 


0=6.4 



(13) 



after having checked that the signal does not change for larger t x , except that the plateaus 
become slightly shorter. To account for the multiplicative renormalisation we proceed as 
follows. We employ the method described in detail in ref. ^U], to compute the renormalisa- 
tion constants Z(afi) and Z{a\x) in the RI/MOM scheme at about 20 different values of the 
scale afj,. We then convert such renormalised ratios R staind '(ty) and R w ^° snbtr '(t y ) into their 
renormalisation invariant forms by using the perturbative anomalous dimension known at 
NLO accuracy [S], namely 

(GO* = c^r 70/2/3 ° (l + ^^ri/mom) (QM) R1/MOM , (14) 
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Figure 1: Plateaus for the ratios defined in eq. J3J) and \1'J\) which are referred to as the standard 
procedure (right plots) and the one without subtractions (left plots). Plateaus are displayed for 
all 4 quark masses and for all three /3's used in this work. 



where 70 = 4 and (3q = 11 — 2nf/3 are universal and 

17397 - 2070n f + 104n? 
Jhi/mom = 8 log 2 6(33 _ 2wf)2 ■ (15) 

The plateaus used to fit the ratios R stand (t y ) and R w/o suhtI -(t y ) to constants R stan± and 
R w /° subtr - ; respectively, for each value of the renormalisation scale a/i, are 

plateau G { [36, 42] /3=6 . , [43, 52)^6.2, [49, 58] /3=6 . 4 } . (16) 

In fig. 121 we illustrate both ratios computed at 24 different values of the renormalisation 
scale at /3 — 6.2, and then converted to the renormalisation group invariant form. After 
proceeding similarly for the other lattice spacings, we find that for (ro/x) 2 > 40, the ratios 
R sta - nd - and R w /° subtr ' nicely follow the perturbatively established scale dependence ([14)1 . 
After fitting these results to a constant on the interval 40 < (Vo/i) 2 < 100 we obtain the 
results listed in table [TJ To convert from a/i to r /i we used the accurately estimated r /a 
from ref. [13], while in the calculation of the two-loop running coupling, a s (r /x), we used 
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Figure 2: Fit of the ratio R{t y ) [see eqs. i)5IJ2)) / to a constant R = i? rgl in the interval 40 < 
(//ro) 2 < 100. The lower plot corresponds to the case in which the effect of mixing with other 
AS = 2 operators has been subtracted. In the upper plot the Ward identity method (without 
subtractions) has been employed. Illustration is provided with the results collected at (3 = 6.2, 
and for k = 0.1344. 



the quenched value r A^" ' = 0.602(48) 05]. In the same tabled we also give the values 
of the bare quark mass p, computed from the axial Ward identity (JTOj) . and of the quantity 
X{m q ) defined as 

8 Z\ (A (t)Al(0)) 8 p P ml _ 

3 (K° P (t)K°;(0)) 3|(0|^|rr°)| 2 " 1 qh { > 

where mp and fp are the mass and the decay constant of the pseudoscalar meson consisting 
of two degenerate quarks of mass m q . 

Before discussing our results for Bk, two important remarks are in order. In the calcula- 
tion of _R stand - we used the subtraction constants given in ref. ^U] where beside the statistical 
we also quoted the systematic uncertainties which arise from the spread of values obtained 
at various values of the momenta flowing through external legs of the elementary vertices 
(i.e., various ap in the RI/MOM scheme). In this paper the subtraction and renormalisation 
constants computed in ^U] are combined with bare matrix elements. The above mentioned 
systematic uncertainties are accounted for by computing the renormalised and subtracted 
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K 


p x 10 2 


X{m q ) 


_R stand -(a) 


R w /° subtr -(a) 


6.0 


0.1335 


5.997(7) 


0.1256(20) 


0.1356(59) 


0.1152(45) 




0.1338 


4.368(8) 


0.1051(19) 


0.1124(55) 


0.0925(42) 




0.1340 


3.750(8) 


0.0910(18) 


0.0969(52) 


0.0776(39) 




0.1342 


3.129(9) 


0.0765(16) 


0.0812(48) 


0.0629(36) 


6.2 


0.1339 


5.792(7) 


0.1766(53) 


0.188(13) 


0.1732(97) 




0.1344 


4.268(7) 


0.1348(48) 


0.142(12) 


0.1277(86) 




0.1349 


2.748(7) 


0.0891(40) 


0.092(12) 


0.0802(72) 




0.1352 


1.834(8) 


0.0589(35) 


0.060(11) 


0.0503(57) 


6.4 


0.1347 


3.144(2) 


0.1336(45) 


0.1363(93) 


0.1226(79) 




0.1349 


2.540(3) 


0.1078(45) 


0.1087(88) 


0.0954(74) 




0.1351 


1.937(3) 


0.0809(44) 


0.0804(79) 


0.0681(65) 




0.1353 


1.334(3) 


0.0536(41) 


0.0518(65) 


0.0409(49) 



Table 1: The values of the (bare) quark mass, p, obtained by using the axial Ward identity ifTTJj) . of 
the quantity X, defined in eq. ijl 7)) . and of the matrix element R obtained by using both methods. 
R are computed from the fit of the ratios and i)I2j) to a constant on the plateau intervals 
indicated in eq. lfW\) . 

matrix element for each value of ap separately. Second important remark is related to the 
accuracy of the two methods employed in this paper. From table Q we see that the errors 
in i? w /° subtr (a) and in i? stand '(a) are comparable. As mentioned in introduction the com- 
putation of the 4-point correlation function needed for i? W// ° subtr (a) is more demanding so 
that -even though one avoids making subtractions- the statistical errors of the two meth- 
ods are essentially equal. The benefit of the method without subtraction is therefore not 
in improving the statistical quality of the results but rather in preventing the occurence 
of uncontrollable systematic uncertainties that might plague the standard method due to 
delicate cancellations of subtractions. 
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a/r 


Rstand. 
K 


A w/o subtr. 
K 


6.0 


0.1863 


1.119(54) 


1.066(39) 


6.2 


0.1354 


1.074(49) 


1.041(37) 


6.4 


0.1027 


1.058(44) 


1.017(46) 


oo 





0.980(114) 


0.961(103) 



Table 2: Results for the Bk parameter as obtained through the fit for all three values of 
the lattice spacing and by using both strategies (standard and the one without subtractions). The 
values of a/rQ are taken from ref. \15jj . We also show the results of the linear extrapolation in 
lattice spacing to the continuum limit. 



3.2 Bx 

With Wilson fermions, 0(a) lattice artifacts can affect the chiral behaviour of the matrix 
elements relevant to the computation of Bk- A convenient way for a clean extraction of Bk 
has been explained in ref. [HJ and consists in studying the dependence of the ratios R stand - 
and i? w /° subtr ' as functions of X, namely 

R = a + (3 X , (18) 

where the fit parameter (3 is identified as Bk(o), and a is the parameter that measures a 
goodness of the chiral behavior of the ratios R. We find that a for all our lattices is consistent 
with zero. From such fits, at each lattice spacing, we thus obtain Bk{cl), all of which are 
listed in table [21 In the same table we also present the results of the extrapolation to the 
continuum limit. That extrapolation has been made linearly since none of the operators 
used in eqs. ()5I10I12|) has been improved. We see that the results of the two procedures 
lead to a consistent value in the continuum limit. That may also be viewed as an evidence 
that no uncontrollable systematic uncertainty used in the standard method occured while 
performing the delicate subtraction procedure. If we imposed the two methods to produce 
exactly the same result in the continuum limit (similar to what has been done in ref. |18j). 
we would have obtained 

B K = 0.969(67) . (19) 

Our errors, after extrapolating to the continuum limit are quite large anyway and we do 
not attempt to include the quadratic term in the continuum extrapolation. The physical 
volume of all our lattices is about (1.7 fm) 3 . By using the formulae of ref. ^Hl it turns that 
for the pseudoscalar mesons consisting of degenerate quarks and with mass mp > 500 MeV, 
the finite volume effects are negligibly small. In the realistic situation, however, one of the 
valence quarks is the strange one (which we can work with directly on the lattice) and the 
other is <i-quark. That situation would necessitate the chiral extrapolation, which in the 
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Figure 3: Extrapolation to the continuum limit. Empty symbols correspond to the results obtained 
at fixed lattice spacing, whereas the rilled ones are the results of the linear extrapolations. The 
shapes of the symbols correspond to two different strategies to compute Bk, as indicated in the 
legend. 



quenched theory would suffer from the (divergent) quenched chiral logarithms. To assess 
some uncertainty due to the degeneracy we may take the relative difference between the 
chiral logarithmic part known in the degenerate and non-degenerate case in full ChPT. 
With A x = 1 GeV, we obtain that B K for the kaon with non-degenerate quarks would 
be only 2% smaller than the one with degenerate quarks. Finally since our calculations 
are made in the quenched approximation, our result cannot make impact on the world 
average value for B^, which is actually completely dominated by the errors due to the use 
of quenched approximation [20] • 3 It is worth mentioning that the short distance piece in the 
unquenched scenario would lead to Bk larger by only 1% 4- 2% compared to the quenched 
one. Such an estimate arises after replacing rif = by m = 4 in eq. (fUjl and in and 
by using A§1 =4) = 294±^ MeV HQ- 



4 Conclusion 

In this letter we presented the results for the renormalisation group invariant bag parameter, 
B K , computed on the lattice with Wilson quarks. Besides the standard procedure, which 
requires a delicate subtraction of the spurious mixing with other AS = 2, dimension-six, 

3 A complete list of results for Bk by using other quark actions with recently updated references can be 
found in ref. [H]. 
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four-quark operators, we also implemented the method based on the use of a Ward identity 
that allows us to avoid the subtraction procedure altogether. 

Our lattice data are produced in the quenched approximation at three values of the 
lattice spacing. At each lattice spacing we use the non-perturbatively computed renormal- 
isation and subtraction constants, presented in ref. ^H]. The conversion to the standard 
renormalisation invariant form is made after checking that our data follow the renormalisa- 
tion scale dependence described by the RI/MOM anomalous dimension coefficients known 
to two-loops in perturbation theory. After having extrapolated to the continuum limit we 
obtain the physically relevant results quoted in table El and eq. As our final estimate 

we chose to quote the results obtained using the method "without subtractions" , namely 
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Appendix 

In this appendix we show that the last two terms of eq. (jEJ) vanish in the SU(3) symmetric 
limit m u = = m s = m. 

We will use charge conjugation and 75-hermiticity which, on the quark propagators 
Sf(x, y; U) (/ is the flavour and U the background gauge configuration), act in the following 
way: 



where the superscripts T and f indicate respectively the transpose and the hermitian con- 
jugation on color and dirac indices. 

Using these two symmetry properties it is easy to show that the trace of an arbitrary 
number of quark propagators and matrices Ti G {I, 75, 7 M , 7 M 7s, <V = °\w = 

75(7^}, computed on a gauge configuration U c , is the complex conjugate of that computed 
on the gauge configuration U, i.e., 



B K = 0.96(10) . 



(20) 



charge conj. C S f (x, y; U) = 'jo^Sjiy, x; U c )^ 2 lo 
hermitian conj. H Sf(x,y; 



(21) 



Tr [T 1 S 1 (x 1 ,x 2 ; U)T 2 S 2 (x 2 , x 3 ; U) . . . T n S n (x n , x\\ U)] = 
Tr [riS , i(xi,x 2 ; U c )T 2 S 2 (x 2 , x 3 ; U c ) . . . T n S n (x n , Xi] U c )]* 



(22) 
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This means that taking the real part of the trace corresponds to the inclusion of the charge- 
conjugated configuration U c in the gauge average. Since the QCD action is symmetric 
under the charge conjugation, the average over iV con f. — > oo will contain the average over 
the configuration U and its charge-conjugated one U c . 

Another property needed is easily obtained by using hermitian conjugation and reads 

Tr [T 1 S 1 (x 1 , x 2 ; U)r 2 S 2 (x 2 , x 3 ; U) . . . T n S n (x n , x x ; U)] = 
[n^(r 4 )] Tr [S n (x 1} x n ,U c )T n ...S 2 (x 3} x 2 ;U c )T 2 S 1 (x 2 ,x 1 ;U c )T 1 ] (23) 

where £(Ti) = +1 for G {1,75,7/^75} and £(Ti) = -1 for Fi G {j^a^,^^}. 

We now analyze the correlation function of Qi between a scalar and a pseudoscalar 
source (since we work in the SU(3) symmetric limit we will not display the flavour indices): 

^(K° s (x)Q 1 (0)K° P (y)) = (d(x)s(x) s(0) 7A .d(0)s(0) 7M 75d(0) d(y) l5 s(y)) = (24) 

(Tr [S(x, 0; U)^S(0, x; U)]Tr [ W S(0, y; U) l5 S(y, 0; U)} 
+ Tr [S{x, 0; U)^ s S{0, x; U)]Tv [ lfl S(0, y; U) lB S(y, 0; U)} 

- Tr [S{x, 0; U) lfl S(0, y; U) lh S{y, 0; U) W S(0, x; U)) 

- Tr [S(x, 0; U) W S{0, y; U)-y s S{y, 0; C/) 7 ^(0, x; U)]) v , 

where (•••)[/ denotes the average over gauge field configurations. 

Using eq. (j2SJ), we see immediately that the sum of traces in eq. ({231 is equal to the 
same expression computed on U c times ^(7^)^(7^75) = — 1. Thus, including the charge- 
conjugated configurations in the gauge average give identically zero for this correlator. Were 
we not working with degenerate m s and m u = masses, these terms should be exponen- 
tially suppressed with respect to the kaon contribution in the limit of large time distances, 
because they correspond to the propagation of scalar states. This point can be explicitly 
checked by computing (Kg{t x )Qi{0)Kp{t y )) in the same numerical simulation for the other 
correlation functions appearing in eq. (JSJ). 
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